clear
clear matrix
set matsize 11000, permanently
set maxvar 32767 
set matsize 800
********************************************************************************

use people
sort id_prof anno
merge 1:1 id_prof anno using nisi /*nisi.dta comes from the WoK database */
tab _merge
drop _merge
sort id_prof anno
merge 1:1 id_prof anno using hindex /*hindex generated using bibliom.do*/
tab _merge
drop _merge
encode ssd, g(ssdcode)
compress

****************************************************************************************
* imputing production to newass and neword - creation of relevant variables for periods pre-entrance
egen begincareer=min(anno), by(id_prof)
egen endcareer=max(anno), by(id_prof)
egen outass=max(newass), by(id_prof)
egen outord=max(neword), by(id_prof)
egen intass=max(ricass), by(id_prof)
egen intord=max(assord), by(id_prof)

* there are strange cases of interrupted careers (probably leaves for different reasons)
tab id_prof if outass==1&outord==1
sort id_prof anno
br id_prof anno ssd fascia cognome nome newass neword if outass==1&outord==1

egen tmp=max(female), by(id_prof)
drop female
ren tmp female
egen tmp=max(birthyear), by(id_prof)
g tmp1=anno-tmp
replace age=tmp1 if age==.
drop tmp*
egen tmp=max(ssdcode), by(id_prof)
replace ssdcode=tmp if ssdcode==.
gen age2=age^2
xtset id_prof anno

* definining the estimation sample, restricting to newass+ricass in 5 year interval after being selected - similarly for full professors
g afterpromass=anno if ricass==1|newass==1
g afterpromord=anno if assord==1|neword==1
g beforepromass=anno if ricass==1|newass==1
g beforepromord=anno if assord==1|neword==1

g tmpass=anno if ricass==1|newass==1
g tmpord=anno if assord==1|neword==1
replace tmpass=tmpass+4
replace tmpord=tmpord+4
egen tmp2ass=max(tmpass), by(id_prof)
egen tmp2ord=max(tmpord), by(id_prof)
replace tmpass=tmpass-8
replace tmpord=tmpord-8
egen tmp3ass=max(tmpass), by(id_prof)
egen tmp3ord=max(tmpord), by(id_prof)

replace afterpromass=(l.afterpromass+1) if afterpromass==.&(l.afterpromass+1)<=(tmp2ass)
replace afterpromord=(l.afterpromord+1) if afterpromord==.&(l.afterpromord+1)<=(tmp2ord)
replace beforepromass=(f.beforepromass-1) if beforepromass==.&(f.beforepromass-1)>=(tmp3ass)
replace beforepromord=(f.beforepromord-1) if beforepromord==.&(f.beforepromord-1)>=(tmp3ord)
replace beforepromass=(f.beforepromass-1) if beforepromass==.&(f.beforepromass-1)>=(tmp3ass)
replace beforepromord=(f.beforepromord-1) if beforepromord==.&(f.beforepromord-1)>=(tmp3ord)
replace beforepromass=(f.beforepromass-1) if beforepromass==.&(f.beforepromass-1)>=(tmp3ass)
replace beforepromord=(f.beforepromord-1) if beforepromord==.&(f.beforepromord-1)>=(tmp3ord)
replace beforepromass=(f.beforepromass-1) if beforepromass==.&(f.beforepromass-1)>=(tmp3ass)
replace beforepromord=(f.beforepromord-1) if beforepromord==.&(f.beforepromord-1)>=(tmp3ord)

drop tmp*

* descriptive statistics - newass and new ord are clearly more productive in the 5 years after promotions, though they have a lower H index
* since they start from zero

replace output1e=0 if output1e==.
replace output1b=0 if output1b==.

sum output1e rH if intass==1&afterpromass!=.
sum output1e rH if outass==1&afterpromass!=.

sum output1e rH if intord==1&afterpromord!=.
sum output1e rH if outord==1&afterpromord!=.

* using individual fixed effect

xtreg output1e age age2 if (intass==1|outass==1)&(afterpromass!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromass!=.&outass==1, xb
gen foutput1e=tmp2+tmp3
* one should use truncreg is one wants to avoid this
replace foutput1e=0 if foutput1e<0&foutput1e!=.
drop tmp*

xtreg output1b age age2 if (intass==1|outass==1)&(afterpromass!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromass!=.&outass==1, xb
gen foutput1b=tmp2+tmp3
* one should use truncreg is one wants to avoid this
replace foutput1b=0 if foutput1b<0&foutput1b!=.
drop tmp*

xtreg rH age age2 if (intass==1|outass==1)&(afterpromass!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromass!=.&outass==1, xb
gen frHe=tmp2+tmp3
* one should use truncreg is one wants to avoid this
replace frH=0 if frH<0&frH!=.
drop tmp*

xtreg output1e age age2 if (intord==1|outord==1)&(afterpromord!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromord!=.&outord==1, xb
replace foutput1e=tmp2+tmp3 if foutput1e==.&beforepromord!=.&outord==1
* one should use truncreg is one wants to avoid this
replace foutput1e=0 if foutput1e<0&foutput1e!=.
drop tmp*

xtreg output1b age age2 if (intord==1|outord==1)&(afterpromord!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromord!=.&outord==1, xb
replace foutput1b=tmp2+tmp3 if foutput1b==.&beforepromord!=.&outord==1
* one should use truncreg is one wants to avoid this
replace foutput1b=0 if foutput1b<0&foutput1b!=.
drop tmp*

xtreg rH age age2 if (intord==1|outord==1)&(afterpromord!=.), fe
predict tmp1, u
egen tmp2=max(tmp1), by(id_prof)
predict tmp3 if beforepromord!=.&outord==1, xb
replace frH=tmp2+tmp3 if frH==.&beforepromord!=.&outord==1
* one should use truncreg is one wants to avoid this
replace frH=0 if frH<0&frH!=.
drop tmp*

sum output1e output1b rH if intass==1&beforepromass!=.
sum foutput1e foutput1b frH if outass==1&beforepromass!=.

sum output1e output1b rH if intord==1&afterpromord!=.
sum foutput1e foutput1b frH if outord==1&afterpromord!=.

* replacing unobserved vaues with predicted ones for newass and neword for 5 years

replace output1e=foutput1e if outass==1&beforepromass!=.&newass!=1
replace output1e=foutput1e if outord==1&beforepromord!=.&neword!=1
replace output1b=foutput1b if outass==1&beforepromass!=.&newass!=1
replace output1b=foutput1b if outord==1&beforepromord!=.&neword!=1

************timing by periods**********************************************************
gen period=1 if anno>=1990&anno<1995
replace period=2 if anno>=1995&anno<1999
replace period=3 if anno>=1999&anno<2003
replace period=4 if anno>=2003&anno<2007
replace period=5 if anno>=2007&anno<2012
label define period 1 "1990-94" 2 "1995-98" 3 "1999-02" 4 "2003-06" 5 "2007-11"

********************* creating the cumulative value of output1 **************************
* we need to drop missing observations for people not yeat entered in academia
tab missing
replace missing=0 if beforepromass!=.
replace missing=0 if beforepromord!=.
tab missing
drop if missing==1

bysort id_prof (anno): gen cumoutput1b=sum(output1b)
bysort id_prof (anno): gen cumoutput1e=sum(output1e) 
 
corr cumoutput1e rH cumoutput1b
bysort subject: corr cumoutput1e rH cumoutput1b if subject!=.
bysort subject: corr cumoutput1e rH cumoutput1b if fascia==1&subject!=.
bysort subject: corr cumoutput1e rH cumoutput1b if fascia==2&subject!=.
bysort subject: corr cumoutput1e rH cumoutput1b if fascia==3&subject!=.

***********************creating two alternative source of ranking*************************
factor cumoutput1e rH , pcf
predict productE

factor cumoutput1b rH , pcf
predict productB
sum productE productB
corr productE productB
save tmp, replace
************************ construction of M index for associate competition - each period is one competition

preserve
egen winass=sum(ricass+newass) if fascia==1|ricass==1|newass==1, by(ssd period)
drop if winass==0|winass==.
egen concorsoass=group(ssd period) if fascia==1|ricass==1|newass==1
sum concorsoass

**** productivity is the final year of the period
collapse subject (max) mricass=ricass mnewass=newass mproductE=productE  mproductB=productB concorsoass area, by(id_prof ssd period)
ren mricass ricass
ren mnewass newass
ren mproductE productE
ren mproductB productB

**** number of participants, positions and winners by concorsoass
g tmp=1
egen N=sum(tmp) if concorsoass!=., by(concorsoass)
egen K=sum(ricass+newass) if concorsoass!=., by(concorsoass)
egen tmp1=sum(newass) if concorsoass!=., by(concorsoass)
gen extass=tmp1/K if concorsoass!=.
drop tmp1

**** number of participants and positions by concorsoord - missing values in the index when K=N

forvalues x = 1/50 {
set seed `x'
gen double tproductE = productE+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductE ) if concorsoass!=., field by(concorsoass)
egen D=sum(R*(ricass+newass)) if concorsoass!=., by(concorsoass)
gen MassE`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoass!=.
drop R D tproductE
}

sum MassE*
corr MassE*
egen mMassE=rowmean(MassE1-MassE50)
egen sdMassE=rowsd(MassE1-MassE50)

forvalues x = 1/50 {
set seed `x'
gen double tproductB = productB+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductB ) if concorsoass!=., field by(concorsoass)
egen D=sum(R*(ricass+newass)) if concorsoass!=., by(concorsoass)
gen MassB`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoass!=.
drop R D tproductB
}

sum MassB*
corr MassB*
egen mMassB=rowmean(MassB1-MassB50)
egen sdMassB=rowsd(MassB1-MassB50)
ren N Nass
ren K Kass

collapse mMassE sdMassE mMassB sdMassB Nass Kass concorsoass extass subject area, by(ssd period)
label value subject subject
label value period period
sort ssd period
save Mass2, replace
restore

set more off 
*********************************replication in order to obtain an initial period including a longer period 1991-98**************************
preserve
recode period 1 2 = 0
label define period 0 "1990-98", modify
drop if period!=0

egen winass=sum(ricass+newass) if fascia==1|ricass==1|newass==1, by(ssd period)
drop if winass==0|winass==.
egen concorsoass=group(ssd period) if fascia==1|ricass==1|newass==1
sum concorsoass

**** productivity is the final year of the period
collapse subject (max) mricass=ricass mnewass=newass mproductE=productE  mproductB=productB concorsoass area, by(id_prof ssd period)
ren mricass ricass
ren mnewass newass
ren mproductE productE
ren mproductB productB

**** number of participants, positions and winners by concorsoass
g tmp=1
egen N=sum(tmp) if concorsoass!=., by(concorsoass)
egen K=sum(ricass+newass) if concorsoass!=., by(concorsoass)
egen tmp1=sum(newass) if concorsoass!=., by(concorsoass)
gen extass=tmp1/K if concorsoass!=.
drop tmp1
**** number of participants and positions by concorsoord - missing values in the index when K=N

forvalues x = 1/50 {
set seed `x'
gen double tproductE = productE+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductE ) if concorsoass!=., field by(concorsoass)
egen D=sum(R*(ricass+newass)) if concorsoass!=., by(concorsoass)
gen MassE`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoass!=.
drop R D tproductE
}

sum MassE*
corr MassE*
egen mMassE=rowmean(MassE1-MassE50)
egen sdMassE=rowsd(MassE1-MassE50)

forvalues x = 1/50 {
set seed `x'
gen double tproductB = productB+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductB ) if concorsoass!=., field by(concorsoass)
egen D=sum(R*(ricass+newass)) if concorsoass!=., by(concorsoass)
gen MassB`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoass!=.
drop R D tproductB
}

sum MassB*
corr MassB*
egen mMassB=rowmean(MassB1-MassB50)
egen sdMassB=rowsd(MassB1-MassB50)
ren N Nass
ren K Kass

collapse mMassE sdMassE mMassB sdMassB Nass Kass extass concorsoass subject area, by(ssd period)
label value subject subject
label value period period
sort ssd period
save Mass2append, replace
use Mass2, clear
append using Mass2append
label define period 0 "1990-98", modify
save Mass2extended, replace
tab period
restore


***********************************************************************************************************************************************

preserve
egen winord=sum(assord+neword) if fascia==2|assord==1|neword==1, by(ssd period)
drop if winord==0|winord==.
egen concorsoord=group(ssd period) if fascia==2|assord==1|neword==1
sum concorsoord

**** productivity is the final year of the period
collapse subject (max) massord=assord mneword=neword mproductE=productE  mproductB=productB concorsoord area, by(id_prof ssd period)
ren massord assord
ren mneword neword
ren mproductE productE
ren mproductB productB


**** number of participants and positions by concorsoord - missing values in the index when K=N
g tmp=1
egen N=sum(tmp) if concorsoord!=., by(concorsoord)
egen K=sum(assord+neword) if concorsoord!=., by(concorsoord)
egen tmp1=sum(neword) if concorsoord!=., by(concorsoord)
gen extord=tmp1/K if concorsoord!=.
drop tmp1

forvalues x = 1/50 {
set seed `x'
gen double tproductE = productE + (-1 + 2*runiform() )/10000
egen R=rank(tproductE ) if concorsoord!=., field by(concorsoord)
egen D=sum(R*(assord+neword)) if concorsoord!=., by(concorsoord)
gen MordE`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoord!=.
drop R D tproductE
}

sum MordE*
corr MordE*
egen mMordE=rowmean(MordE1-MordE50)
egen sdMordE=rowsd(MordE1-MordE50)

forvalues x = 1/50 {
set seed `x'
gen double tproductB = productB+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductB ) if concorsoord!=., field by(concorsoord)
egen D=sum(R*(assord+neword)) if concorsoord!=., by(concorsoord)
gen MordB`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoord!=.
drop R D tproductB
}

sum MordB*
corr MordB*
egen mMordB=rowmean(MordB1-MordB50)
egen sdMordB=rowsd(MordB1-MordB50)
ren N Nord
ren K Kord

collapse mMordE sdMordE mMordB sdMordB Nord Kord concorsoord extord subject area, by(ssd period)
label value subject subject
label value period period
sort ssd period
save Mord2, replace
restore

*********************************replication in order to obtain an initial period including a longer period 1991-98**************************
preserve
recode period 1 2 = 0
label define period 0 "1990-98", modify
drop if period!=0

egen winord=sum(assord+neword) if fascia==2|assord==1|neword==1, by(ssd period)
drop if winord==0|winord==.
egen concorsoord=group(ssd period) if fascia==2|assord==1|neword==1
sum concorsoord

**** productivity is the final year of the period
collapse subject (max) massord=assord mneword=neword mproductE=productE  mproductB=productB concorsoord area, by(id_prof ssd period)
ren massord assord
ren mneword neword
ren mproductE productE
ren mproductB productB


**** number of participants and positions by concorsoord - missing values in the index when K=N
g tmp=1
egen N=sum(tmp) if concorsoord!=., by(concorsoord)
egen K=sum(assord+neword) if concorsoord!=., by(concorsoord)
egen tmp1=sum(neword) if concorsoord!=., by(concorsoord)
gen extord=tmp1/K if concorsoord!=.
drop tmp1

forvalues x = 1/50 {
set seed `x'
gen double tproductE = productE + (-1 + 2*runiform() )/10000
egen R=rank(tproductE ) if concorsoord!=., field by(concorsoord)
egen D=sum(R*(assord+neword)) if concorsoord!=., by(concorsoord)
gen MordE`x'=2*((((K*(2*N-K+1))/2)-D)/(K*(N-K)))-1 if concorsoord!=.
drop R D tproductE
}

sum MordE*
corr MordE*
egen mMordE=rowmean(MordE1-MordE50)
egen sdMordE=rowsd(MordE1-MordE50)

forvalues x = 1/50 {
set seed `x'
gen double tproductB = productB+ (-1 + 2*runiform() )/10000000
egen R=rank(tproductB ) if concorsoord!=., field by(concorsoord)
egen D=sum(R*(assord+neword)) if concorsoord!=., by(concorsoord)
gen MordB`x'=(((K*(2*N-K+1))/2)-D)/(K*(N-K)) if concorsoord!=.
drop R D tproductB
}

sum MordB*
corr MordB*
egen mMordB=rowmean(MordB1-MordB50)
egen sdMordB=rowsd(MordB1-MordB50)
ren N Nord
ren K Kord

collapse mMordE sdMordE mMordB sdMordB Nord Kord concorsoord extord subject area, by(ssd period)
label value subject subject
label value period period
sort ssd period
save Mord2append, replace
use Mord2, clear
append using Mord2append
label define period 0 "1990-98", modify
save Mord2extended, replace
tab period
restore

* descriptive statistics of the Mindex - unit of analysis are the ssd and period - old results
use Mass2, clear
g KsuN=Kass/Nass
table period subject, c(m mMassE sd mMassE m mMassB sd mMassB m KsuN) row col format(%9.2f)

use Mord2, clear
g KsuN=Kord/Nord
table period subject, c(m mMordE sd mMordE m mMordB sd mMordB m KsuN) row col format(%9.2f)

* descriptive statistics of the Mindex - unit of analysis are the ssd and period - new results

use Mass2extended, clear
gen biblio=(area==1|area==2|area==3|area==4|area==5|area==6|area==7|area==9|ssd=="ICAR/01"|ssd=="ICAR/02"|ssd=="ICAR/03"|ssd=="ICAR/04"|ssd=="ICAR/05"| ///
ssd=="ICAR/06"|ssd=="ICAR/07"|ssd=="ICAR/08"|ssd=="ICAR/09"|ssd=="ICAR/22"| ///
ssd=="M-PSI/01"|ssd=="M-PSI/02"|ssd=="M-PSI/03"|ssd=="M-PSI/04"|ssd=="M-PSI/05"|ssd=="M-PSI/06"|ssd=="M-PSI/07"|ssd=="M-PSI/08")

label define areacun 1 "mathematics and computer science" 2 "physics" 3 "chemistry" 4 "earth science" 5 "biology" ///
6 "medicine" 7 "agriculture and veterinary sciences" 8 "engineering and architecture(biblio)" 15 "architecture (not biblio)" 9 "industrial engineering and ICT" ///
10 "humanities" 11 "history, philosophy" 16 "psichology" 12 "law" 13 "economics and statistics" 14 "sociology and political science" 99 "missing"
label value areacun areacun

g KsuN=Kass/Nass

table period subject if period>0, c(m mMassE sd mMassE m mMassB sd mMassB m extass) row col format(%9.3f)
table area period if period>0&area<15, c(m mMassE sd mMassE m mMassB sd mMassB m extass) row col format(%9.3f)
table period biblio if period>0, c(m mMassE sd mMassE m mMassB sd mMassB m extass) row col format(%9.3f)

gen areabis=areacun
replace areabis=15 if ssd=="ICAR/10"|ssd=="ICAR/11"|ssd=="ICAR/12"|ssd=="ICAR/13"|ssd=="ICAR/14"|ssd=="ICAR/15"|ssd=="ICAR/16"|ssd=="ICAR/17"|ssd=="ICAR/18"| ///
ssd=="ICAR/19"|ssd=="ICAR/20"|ssd=="ICAR/21"
replace areabis=16 if ssd=="M-PSI/01"|ssd=="M-PSI/02"|ssd=="M-PSI/03"|ssd=="M-PSI/04"|ssd=="M-PSI/05"|ssd=="M-PSI/06"|ssd=="M-PSI/07"|ssd=="M-PSI/08"
label value areabis areacun

**** main table 
table areabis period if period>0&areabis<20, c(m mMassE sd mMassE n mMassE) row col format(%9.3f)

drop if period==0
twoway (kdensity mMassE if areabis==1) (kdensity mMassE if areabis==2) (kdensity mMassE if areabis==3) (kdensity mMassE if areabis==4) ///
(kdensity mMassE if areabis==5) (kdensity mMassE if areabis==6) (kdensity mMassE if areabis==7) (kdensity mMassE if areabis==8) /// 
(kdensity mMassE if areabis==9) (kdensity mMassE if areabis==16), title("bibliometric") ytitle("kdensity") saving(tmp1, replace)

twoway (kdensity mMassE if areabis==15) (kdensity mMassE if areabis==10) (kdensity mMassE if areabis==11) (kdensity mMassE if areabis==12) ///
(kdensity mMassE if areabis==13) (kdensity mMassE if areabis==14), title("not bibliometric") ytitle("kdensity") saving(tmp2, replace)
graph combine tmp1.gph tmp2.gph, ti(M index by research field - associate professors)

twoway (kdensity mMassE if biblio==0) (kdensity mMassE if biblio==1)

twoway (kdensity mMassE if areabis==1) (kdensity mMassE if areabis==2) (kdensity mMassE if areabis==3) (kdensity mMassE if areabis==4) ///
(kdensity mMassE if areabis==5) (kdensity mMassE if areabis==6) (kdensity mMassE if areabis==7) (kdensity mMassE if areabis==8) /// 
(kdensity mMassE if areabis==9) (kdensity mMassE if areabis==16)(kdensity mMassE if areabis==15) (kdensity mMassE if areabis==10) (kdensity mMassE if areabis==11) (kdensity mMassE if areabis==12) ///
(kdensity mMassE if areabis==13) (kdensity mMassE if areabis==14), legend(off)

*** full professors

use Mord2extended, clear
gen biblio=(area==1|area==2|area==3|area==4|area==5|area==6|area==7|area==9|ssd=="ICAR/01"|ssd=="ICAR/02"|ssd=="ICAR/03"|ssd=="ICAR/04"|ssd=="ICAR/05"| ///
ssd=="ICAR/06"|ssd=="ICAR/07"|ssd=="ICAR/08"|ssd=="ICAR/09"|ssd=="ICAR/22"| ///
ssd=="M-PSI/01"|ssd=="M-PSI/02"|ssd=="M-PSI/03"|ssd=="M-PSI/04"|ssd=="M-PSI/05"|ssd=="M-PSI/06"|ssd=="M-PSI/07"|ssd=="M-PSI/08")

label define areacun 1 "mathematics and computer science" 2 "physics" 3 "chemistry" 4 "earth science" 5 "biology" ///
6 "medicine" 7 "agriculture and veterinary sciences" 8 "engineering and architecture(biblio)" 15 "architecture (not biblio)" 9 "industrial engineering and ICT" ///
10 "humanities" 11 "history, philosophy" 16 "psichology" 12 "law" 13 "economics and statistics" 14 "sociology and political science" 99 "missing"

g KsuN=Kord/Nord
table period subject if period!=1&period!=2, c(m mMordE sd mMordE m mMordB sd mMordB m extord) row col format(%9.2f)
table area period if period!=1&period!=2&area!=99, c(m mMordE sd mMordE m mMordB sd mMordB m extord) row col format(%9.2f)
table period biblio if period!=1&period!=2, c(m mMordE sd mMordE m mMordB sd mMordB m extord) row col format(%9.2f)

gen areabis=areacun
replace areabis=15 if ssd=="ICAR/10"|ssd=="ICAR/11"|ssd=="ICAR/12"|ssd=="ICAR/13"|ssd=="ICAR/14"|ssd=="ICAR/15"|ssd=="ICAR/16"|ssd=="ICAR/17"|ssd=="ICAR/18"| ///
ssd=="ICAR/19"|ssd=="ICAR/20"|ssd=="ICAR/21"
replace areabis=16 if ssd=="M-PSI/01"|ssd=="M-PSI/02"|ssd=="M-PSI/03"|ssd=="M-PSI/04"|ssd=="M-PSI/05"|ssd=="M-PSI/06"|ssd=="M-PSI/07"|ssd=="M-PSI/08"
label value areabis areacun

**** main table 
table areabis period if (period!=1|period!=2)&areabis<20, c(m mMordE sd mMordE n mMordE) row col format(%9.3f)

drop if period==1|period==2
twoway (kdensity mMordE if areabis==1) (kdensity mMordE if areabis==2) (kdensity mMordE if areabis==3) (kdensity mMordE if areabis==4) ///
(kdensity mMordE if areabis==5) (kdensity mMordE if areabis==6) (kdensity mMordE if areabis==7) (kdensity mMordE if areabis==8) /// 
(kdensity mMordE if areabis==9) (kdensity mMordE if areabis==16), title("bibliometric") ytitle("kdensity") saving(tmp1, replace)

twoway (kdensity mMordE if areabis==15) (kdensity mMordE if areabis==10) (kdensity mMordE if areabis==11) (kdensity mMordE if areabis==12) ///
(kdensity mMordE if areabis==13) (kdensity mMordE if areabis==14), title("not bibliometric") ytitle("kdensity") saving(tmp2, replace)
graph combine tmp1.gph tmp2.gph, ti(M index by research field - full professors)

twoway (kdensity mMordE if biblio==0) (kdensity mMordE if biblio==1)

twoway (kdensity mMordE if areabis==1) (kdensity mMordE if areabis==2) (kdensity mMordE if areabis==3) (kdensity mMordE if areabis==4) ///
(kdensity mMordE if areabis==5) (kdensity mMordE if areabis==6) (kdensity mMordE if areabis==7) (kdensity mMordE if areabis==8) /// 
(kdensity mMordE if areabis==9) (kdensity mMordE if areabis==16)(kdensity mMordE if areabis==15) (kdensity mMordE if areabis==10) (kdensity mMordE if areabis==11) (kdensity mMordE if areabis==12) ///
(kdensity mMordE if areabis==13) (kdensity mMordE if areabis==14), legend(off)




